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Abstract 


This  paper  combines  recent  developments  in  the  area  of  generation  of  dependent  random 
variables  with  the  advantages  of  the  use  of  common  and  antithetic  random  numbers.  This 
combination  yields  new  efficient  methods  for  estimating  complicated  stochastic  quantities  by 
simulation.  Some  theoretical  and  practical  aspects  of  use  of  antithetic  and  common  random 
numbers  for  variance  reduction  while  using  the  total  hazard  construction  are  given.  A  proof 
of  their  optimality  in  estimating  the  expected  value  of  the  response  sum  or  the  response 
difference  of  functions  of  vector  arguments  with  dependent  components  is  presented.  Some 
numerical  examples  illustrate  the  theory. 

Key  words  and  phrases.  Antithetic  variates,  multivariate  dependence,  Monte  Carlo  methods, 
simulation,  generation  of  dependent  variables,  total  hazard  construction,  coherent  life  functions, 
reliability  theory. 
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Abstract 


1.  Introduction 


Methods  of  variance  reduction  in  simulation  analysis  are  useful  for  the  purpose  of  improving 
accuracy  of  estimates  of  output  parameters  of  complex  stochastic  systems.  The  development 
of  such  methods  is  important  because  they  yield  estimators  which  are  superior  in  accuracy 
to  the  crude  Monte  Carlo  estimator. 

One  of  the  most  successful  variance  reduction  techniques  is  the  use  of  common  and  antithetic 
random  numbers  (see  the  bibliography  of  Rubinstein,  Samorodnitsky  and  Shaked  (1985)). 

This  paper  presents  new  ways  of  application  of  common  and  antithetic  random  variables.  It 
combines  recent  developments  in  the  area  of  generation  of  dependent  random  variables  with 
the  advantages  of  the  use  of  common  and  antithetic  random  variates.  This  combination  yields 
new  efficient  methods  for  estimating  complicated  probabilistic  quantities  by  simulation  (see 
Section  4). 

More  explicitly,  in  the  present  paper  these  methods  are  used  in  order  to  minimize  the  variance 
of  g(X)  -  h(Y)  where  g  and  h  are  real  measurable  functions  which  are  monotonic  in  the  same 
(or  the  opposite)  direction  and  the  joint  cdf  (=cumulative  distribution  function)  of  X  and  Y 
is  restricted  to  belonging  to  some  set  ©  of  cdf’s.  To  see  an  instance  of  such  a  set  0,  suppose, 
for  example,  that  X  =  (ATi, ..., X„)  and  Y=  (yj,...,  Yn)  are  both  n-dimensional  random  vectors 
with  respective  cdf’s  Fx  and  Fy ■  Then  ©  can  be  taken  to  be  the  set  of  all  cdf’s  on  R2n  with 
marginals  Fx  and  Fy  or  (see,  e.g.,  Rubinstein,  Samorodnitsky  and  Shaked  (1985))  the  set  of 
all  cdf’s  on  R2n  with  marginals  Fx  and  Fy  such  that  dependence  is  permitted  only  between 
like  components  of  X  and  Y  (components  with  the  same  indices). 
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Usually  X  and  Y  are  generated  by  constructing  them  from  independent  uniform  random 
variables.  In  this  respect,  three  alternative  constructions  have  been  used.  They  are  (i)  standard 
construction  (e.g.,  Arjas  and  Lehtonen  (1978)),  (ii)  dynamic  construction  (e.g.,  Shaked  and 
Shanthikumar  (1985b)),  and  (iii)  total  hazard  construction  (e.g.,  Norros  (1984),  Shaked  and 
Shanthikumar  (1985a)).  Each  of  these  constructions  determine  vector  functions  Ty  and  Ty 
depending,  respectively,  on  the  cdf’s  of  X  and  Y  such  that  xiT^U)  and  Y=Ty(V)  where  U 
and  V  are  vectors  of  independent  uniform  [0,1]  random  variables  (here  d  denotes  equality 
in  law). 

When  X  and  Y  are  both  n-dimensional  (as  will  be  assumed  throughout  the  paper),  then  the 
standard  and  the  total  hazard  constructions  (to  be  described  below  in  detail)  require  U  and 
V  to  be  n-dimensional.  The  dynamic  construction  requires  U  and  V  to  be  of  (the  same) 
dimension  n  >  n. 

We  say  that  VCRN  (vector  of  common  random  numbers)  is  used  if  V—  U.  We  say  that 
VARN  (vector  of  antithetic  random  numbers)  is  used  if  F=»l-t/ where  1  s  (1, ...,  1)'. 

Let  g  and  h  be  real  measurable  functions  and  suppose  we  are  interested  in  reducing  or 
minimizing  Var(g(X)  —  h(Y))  subject  to  Fyy  €  ©  where  ©  is  some  set  of  2n-dimensional  cdf’s. 
Assume  we  construct  X  and  Y  by  7y((/)  and  Ty(V)  where  U  and  V  are  not  necessarily 
independent  vectors  of  independent  uniform  [0,1]  random  variables.  Then 

(1.1)  Var(g(X)  -  KY))  -  Var{g  *  {U)  —h*  (  F)) 

where 

(1.2)  g*(«)=g(7*(«)), 

(1.3)  h  *  (y)  =  h(Ty(y)). 

i 
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Representation  (1.1)  together  with  the  following  result  (Lemma  1.1  below)  of  Rubinstein  and 
Samorodnitsky  (1982)  [see  also  Rubinstein,  Samorodnitsky  and  Shaked  (1985)]  will  be  used 
to  find  optimal  variance  reduction  techniques  in  Sections  2  and  3.  First  we  need  to  define 
the  following: 

Let  $  be  the  set  of  2n-dimensional  cdf’s  corresponding  to  random  vectors  (U,  V)  of  uniform 
[0,1]  random  variables  with  n-dimensional  marginals  and 

Fy(v)  =  n?=iv(,0  <  Ui  <  1,  0  <  v,  <  1,  i  =  l,...,n,  with  dependence  permitted  only  between  like 
components  of  u  and  y. 

Lemma  1.1.  Let  g  and  h  be  real  measureable  functions  on  [0,1]"  (a)  If  g  and  h  are 

monotonic  in  the  same  direction  with  respect  to  the  i-th  component,  /=  then 


(1.4)  min  Var(g{U)  -h{V))  =  Var{g(U) -h  {0). 

Fuve  4> 

(b)  If  g  and  h  are  monotonic  in  the  opposite  direction  with  respect  to  the  ;-th  component, 
i=  1 then 

(1.5)  mm  Var(g{U)  -h{V))  =  Var{g{U)-h{\  -  U)). 

Fuv  e  <P 

That  is,  use  of  VCRN  is  optimal  for  problem  (1.4)  and  use  of  VARN  is  optimal  for  problem 

(1.5) . 

In  the  following  sections  we  will  use  (1.1)  and  Lemma  1.1  as  follows.  For  a  given  construction 
Tx(lD  aRd  IyUD  we  will  find  conditions  under  which  g*  of  (1.2)  and  h*  of  (1.3)  are 
monotone  in  the  same  [respectively,  opposite]  direction  and  then  from  Lemma  1.1  it  will 
follow  that  use  of  VCRN  [respectively,  VARN]  is  optimal. 
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In  Section  2  we  deal  with  variance  reduction  while  using  the  standard  construction.  We 
briefly  describe  the  main  results  of  Rubinstein,  Samorodnitsky  and  Shaked  (1985)  and  their 
underlying  ideas,  and  in  doing  this  we  motivate  the  developments  of  the  following  section. 
Section  3  deals  with  variance  reduction  while  using  the  total  hazard  construction.  We  follow 
the  lines  of  Section  2,  thus  obtaining  new  results  which  are  the  main  contributions  of  this 
paper.  Finally,  in  Section  4,  we  present  some  simulation  results  involving  dependent  stochastic 
variables. 
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2.  The  standard  construction  and  variance 
reduction 


First  the  standard  construction  is  described  in  detail.  The  following  notation  is  used.  For 
any  random  vector  Z  =  (Zi,...,Z„),  let 

FZ,UX^,J^  =  -  Zi  I  ^1  =  =  2/-l  } 


and  let 


n. 


According  to  the  standard  construction,  X=  (X\,...,Xn)  and  Y=  ( Kj , ....  Yn)  are  generated  as 
follows: 


(2.1)  Xx^F^(Ux)  =  ax(Ux), 


(2.2)  X(  =  FX^  |  . S’.  jCt;, . U,_X)^U^  =  aSU l»-»  ^).'  =  2,...,n, 


(2.3)  Fi^V,)^?,^), 
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(2.4) 


‘F~}~ 


y,l6,(^) . . K,_,)(^  “  b^V i— - 


where  t/  and  Jf  are  vectors  of  independent  uniform  [0,1]  random  variables.  If  we  define 
Oj(  U)  =  Uj)  and  bj{V)  s  Vj)  then  (2.1)  -  (2.4)  can  be  rewritten  as 


(2.5)  Xt  -<%(£#,/ =l,...,/i, 

(2.6)  1^  =  6,  (20, 


Thus,  the  transformation  7^  and  Ty,  mentioned  in  Section  1,  which  corresponds  to  the 
standard  construction  is  described  in  (2.5)-(2.S). 


A  random  vector  Z  =  (Z\, ....  Z„)  is  said  to  be  CIS  (conditionally  increasing  in  sequence)  if 
P(Zj>2j  \Z\  =2i,...,Z/_i  =  i )  is  nondecreasing  in  z\,  ...,z,_i  for  every  z,,  i—  2,...,n.  Rubinstein, 
Samorodnitsky  and  Shaked  (1985)  proved: 


Lemma  2.1.  If  X  and  Y  are  CIS  then  o,(m)  and  bj{v)  [of  (2.5)-(2.6)]  are  nondecreasing 
functions  of  u  and  v, /  =  l,...,n. 

That  is,  if  X  and  Y  are  CIS  then  the  transformations  Tx  and  Ty  defined  by  (2.5)-(2.6)  are 
monotone  nondecreasing. 

Let  Z  be  the  set  of  all  cdf’s  on  F2n  with  marginals  Fx  and  Fy  such  that  dependence  is 
permitted  only  between  like  components  of  U  and  V  of  (2.5)  and  (2.6).  Consider  the  problem 

(2.7)  minimize  Var(g(X)  -  li{Y))  subject  to  FXy€.  Z. 

From  Lemmas  1.1  and  2.1  one  easily  obtains  (see,  e.g.,  Rubinstein,  Samorodnitsky  and  Shaked 
(1985)): 
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Theorem  2.1.  If  X  and  Y  are  CIS  and  g  and  h  are  monotone  in  the  same  (the  opposite) 
direction  then  use  of  VARN  (VCRN)  yields  an  optimal  solution  for  problem  (2.7). 

Note  that  the  cdf  Fxy »  which  corresponds  to  the  choice  of  independent  U  and  V  in  (2.5)-(2.6), 
belongs  to  2.  Thus,  by  Theorem  2.1,  the  use  of  VCRN  (VARN)  is  preferable  to  the  use  of 
independent  t/and  V.  That  is,  use  of  VCRN  (VARN)  is  superior  to  the  crude  Monte  Carlo 
method.  Furthermore,  use  of  VCRN  or  VARN  requires  generation  of  only  half  as  many 
uniform  [0,1]  random  variables. 
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3.  The  total  hazard  construction  and 
variance  reduction. 


The  total  hazard  construction  can  be  used  to  generate  absolutely  continuous  nonnegative 
random  vectors  (see,  e.g.,  Norros  (1984)  and  Shaked  and  Shanthikumar  (1985a)).  It  is 
particularly  useful  when  the  numerical  inversion  of  the  conditional  distributions,  described  in 
(2.1)-(2.4),  is  too  involved  and  time  consuming,  and  on  the  other  side,  inversion  of  the 
multivariate  hazard  functions,  to  be  defined  below,  is  simple.  An  example  of  such  an  instance 
is  given  in  Section  4. 

Let  X—  (A'j,  be  an  absolutely  continuous  nonnegative  random  vector  to  be  generated. 

We  will  use  the  following  notation. 

prr  J=  c  {1, let  xj  denote  (xji,...,xJk).  If  J  -  then  xj  denotes 

(*,-,,  ...,xin  k).  Let  1_  =  ( 1 , ....  1  )#-  The  length  of  i  will  vary  from  one  formula  to  another,  but 
it  will  be  always  possible  to  determine  it  from  the  expression  in  which  1  appears. 


For  /c  and  ieJ  let  A,(x  |  Xj  =  xj,  Xj  >  x_l)  denote  the  conditional  hazard  rate  of  Xj 


at  time  x  given  that  Xj-xj  and  that  Xj  >  xi,  where  x  >  V  x. 

J  j  e  J 

V  xj  =  0.  Formally,  for  i  e  J, 
jeJ 


ma x{xj.j  €  A-  If  J  =  4>  then 
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(3.1)  \i(x\XJ=xj,XJ>x\) 


s,£To  izp{x<xi^x+^zj=*j'Zj>xA'x*jVjxj 


(J  may  be  empty.)  The  absolute  continuity  of  X  ensures  that  limit  exists.  To  save  space 

we  sometimes  supress  the  condition  Xj>x i  and  just  write  A,(x  |  Xj  —  xj,  •  )  but  the  reader 

should  keep  in  mind  that  means  Xj>x  1  with  x  being  the  same  as  the  first  argument  of 

Xj.  Note  that  A,(x  \Xj  —  xj,m)  is  well  defined  for  all  x  >  V  x:. 

jeJJ 

For  i  eJ  the  total  hazard  accumulated  by  Xj  during  the  time  interval  V  x,-,  V  x,-  +  x  ,  x  >  0, 

Ly  e  J  j  e  J  J 

is  defined  by 


A/0*  I  Xj=  Xj)  =  fJ,JJ  Xjiu]  Xj=  xj,.)du,  x  >0,i  e  J 

J  V  r 


V*l 

JtjJ 


When  J  =  <£,  A,  (x  |  Xj  —  xj)  will  be  simply  denoted  by  A,(x). 

We  will  introduce  now  a  notation  for  the  total  hazard  accumulated  by  Xj  by  time  x.  Fix 
x  >  0  and  suppose  that  it  is  given  that  Ay,, ...,  Xjkl(k  >  1)  failed  at  times  Xjr...,Xjk{,  respectively 
(xy,  < ...  <xJk_l  <x)  and  that  all  other  Ays  are  alive  at  time  x.  For  //{y'i, i }  denote 


*-i 


+  SA^ 'XA  =  =XJ,-? 


1-2 


+A/(*  "*/*-.  I  Xji==xji'  -'XX-l=xX-? 
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Also  denote  (corresponding  to  the  case  k=  1) 

*,(*)  =  A,-(x),jr  >  0. 

The  total  hazard  accumulated  by  Xj  by  the  time  it  f ailed ,  given  that  Xj  was  the  /c-th  to  fail  and 
that  Xjl,...,Xjk_i  failed  before  Xh  is  ^ i\ ju...jk_s(X>  I 

Define  the  inverse  functions 

A"1  (0  =  inf >  0:A,(x)  >  r],  i  =  l,...,n,/  >  0, 

and  for  nonempty  Jc  Jl,. x/>01.  and  i  c  7, 

A”jy(/ 1  xj)  =  in/[x  >  0:A,(x  \  Xj  =  xj)  >  t},t  >0  . 

Let  U\,...,Un  be  independent  uniform  [0,1]  random  variables.  They  generate  independent 
standard  (i.e.,  mean  one)  exponential  random  variables  as  follows: 

(3.3)  Ej  =  —  log(l  —  Uj),  i  =  1, ...» n. 

The  total  hazard  construction  consists  of  transforming  E\,...,En  into  Xi,...,X„  as  follows: 

Step  I.  Let  j\  be  the  (random)  index  (which,  by  absolute  continuity,  is  unique  with  probability 
one)  such  that  A/i 1  (Ejj)  =  min{A,_1  (£,):/ =  l,...,n}  and  define 

<3  4>  \  -  A,T'<V' 
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Step  k  (k=2,...,n).  Given  that  steps  1,2,. ...k-l  resulted  in  j  let  /—  £/j, --yjk—x }•  Let 

jit  be  the  (random)  index  (which,  by  absolute  continuity,  is  unique  with  probability  one)  such 
that 


Jk  \j\' 


i  (*, 
■Jk-l  Jk- 1 


,x  )|x, 

J 1  Jk-l  Jl 


=  min.  { A, "}[£;- 

I  €.  J 


(X:  \X,% 

-1  Jk-l  Jl 


-4JI \ . 4-,]}  • 


It  is  easy  to  verify,  by  induction,  that  the  arguments  of  Af 1  and  A^1  ,  in  the  above  expression, 
are  nonnegative.  Having  chosen  the  (random)  index  as  described  above,  define  (here 

J=  [h . A-i}) 


(3.5) 


X:  =X: 

Jk  Jk-l 


^jk\ji-jkSxjk-i 


X:,...,X:  )|X, 
Jl  Jk-l  Jl 


,n  . 


Notice  that  (3.4)-(3.5)  describe  explicitly  how  to  obtain  X  from  E-  (£j, ...,  En).  However 
(3.4)-(3.5)  together  with  (3.3)  define  a  transformation  of  C/into  X.  This  is  the  transformation 
I*  mentioned  in  Section  1,  which  corresponds  to  the  total  hazard  construction.  Shaked  and 
Shanthikumar  (1985a)  proved  that  indeed  (3.4)-(3.5)  yields  a  random  vector  with  the  desired 
cdf.  Analogously,  starting  from  independent  uniform  [0,1]  random  variables  V\,...,Vn  one 
can  generate  a  nonnegative  absolutely  continuous  random  vector  (yj,...,  Yn)  with  a  desired 
cdf.  That  is,  the  transformation  Ty,  mentioned  in  Section  1,  which  corresponds  to  the  total 
|  hazard  construction,  is  described  in  (3.3)-(3.5)  with  obvious  modification  of  notation. 

■  In  order  to  give  conditions  under  which  and  Ty  are  nondecreasing,  we  need  the  following 

|  definition.  For  2  <  k  <  n  we  use  the  notation 

I 

i 

I, 
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ck c^k*  xU—tX/c_i) 


sxi k-1  Ui.—  I  *i»  •••»•*&_  l] 

which  describes,  according  to  (3.5),  how  A*  is  determined,  given  that  ji  =  1,  ...,Jk-l  —  k  ~  l.jfc  —  ^ 
and  that  Zy.  =  tp,X\  =x\,  j  =  xp—  j. 

Definition  3.1.  The  absolutely  continuous  nonnegative  random  vector  X  =  (X\,...,Xn)  is  said 
to  have  CDTH  {conditionally  decreasing  total  hazard)  if 

ck+2^k+2’Xlf  •"'xk’ck+0'  —  *ci—  l(*i— 1»*1*  "■>xk'ck+ 1»  ci-2^] 

(the  arguments  of  some  of  the  cm’s  are  omitted)  nonincreases  in 

xke  {xk:xk>xk_i,ck+i>xk;  cl+i>cl>t  =  k  +  l . i-2}  for  all  1<£+1  </'<«, 

ft:+l  £  0,  >  0,  0  <  x\  <  ...  and  if  the  above  condition  holds  for  all  permutations 

of  the  indices  1,2 

Definition  3.1  is  essentially  the  same  as  the  definition  of  supportive  system  of  Norros  (1984). 
An  easy  to  check  condition  which  implies  the  CDTH  property  is  given  in  (3.7)  below. 

From  the  proof  of  Theorem  4.4  of  Shaked  and  Shanthikumar  (1985a)  the  following  result 
follows: 

Lemma  3.1.  If  the  absolutely  continuous  nonnegative  random  vectors  X=  (ATj,  ...,Xn)  and 
Y=>  (yj,...,  Vn)  have  CDTH  then  the  transformations  Ty  and  Ty  defined  by  (3.3)-(3.5)  are 
monotone  nondecrcasing. 
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Let  z!  be  the  set  of  all  cdf’s  on  R^1  with  marginals  Fx  and  Fy  such  that  dependence  is 
permitted  only  between  like  components  of  U  and  V  of  the  transformations  Tx(U)  and  Ty(V) 
defined  by  (3.3)-(3.5).  Consider  the  problem 

(3.6)  minimize  Var(g(X)  —  h(Y))  subject  to  FXy€.  z.'  . 

From  Lemmas  1.1  and  3.1  one  easily  obtains 

Theorem  3.1.  If  X  and  Y  have  CDTH  and  g  and  h  are  monotone  in  the  same  (the  opposite) 
direction  then  use  of  VCRN  (VARN)  yields  an  optimal  solution  for  problem  (3.6). 

Note  that  the  cdf  Fxy,  which  corresponds  to  the  choice  of  independent  U  and  V  in  the 
transformation  defined  by  (3.3)-(3.5),  belongs  to  z Thus,  by  Theorem  3.1,  the  use  of 
VCRN  (or  VARN)  is  preferable  to  the  use  of  crude  Monte  Carlo  method  (i.e.,  independent 
U  and  V). 

Shaked  and  Shanthikumar  (1985a)  showed  that  if,  for  a  random  vector  X,  the  multivariate 
conditional  hazard  rates  (defined  in  (3.1))  nondecrease  as  functions  of  the  number  of  past 
failures  and  nonincrease  as  functions  of  the  failure  times  of  these  past  failures,  then  X  has 
CDTH.  That  is,  X  has  CDTH  if  more  and  earlier  failures  cause  a  higher  risk  for  the  surviving 
components  than  fewer  and  later  failures. 

More  explicitly  Shaked  and  Shanthikumar  (1985a)  showed  that  if  for  disjoint  sets  7,/c  Jl,...,n} 
and  fixed  x/,xj,xj  (such  that  xj  <xj)  and  k  e  TUI  (/  or  J  may  be  empty), 

(3.7)  Xka .  V  Xj)  V{ .  V  xj)  +  u  \Xj  =  xj,  Xj  =  xj, .) 

/  €  /  J  €  J 

>K({  YXj)V(  Kx/)  +  u|X/=x/,.),m>0, 

,  i  e  I  j  £  / J 

\ 
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k 


then  (Xi,...,Xn)  has  CDTH.  Thus,  from  Theorem  3.1  we  get 

Corollary  3.1.  If  the  multivariate  hazard  rates  of  X  satisfy  (3.7)  and  the  multivariate  hazard 
rates  of  Y  satisfy  the  same  condition  with  proper  notational  modifications,  and  if  g  and  h 
are  monotone  in  the  same  (the  opposite)  direction  then  use  of  VCRN  (VARN)  yields  an 
optimal  solution  for  problem  (3.6). 


Norros  (1984)  has  obtained  another  condition  which  implies  that  the  transformation  defined 
by  (3.4)-(3.5)  is  monotone  nondecreasing.  Thus,  if  Norros’  condition  holds  for  X  and  Y  and 
if  g  and  h  are  monotone  in  the  same  (the  opposite)  direction  then  use  of  VCRN  (VARN) 
yields  an  optimal  solution  for  problem  (3.6). 
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4.  Some  applications  and  examples 


In  this  section  we  illustrate  the  use  of  the  results  of  Section  3  for  simulation  purposes. 

Suppose  we  want  to  estimate  0*^  =  E[g(X)]  where  X  has  CDTH  and  g  is  a  real  valued 
measurable  monotonic  function.  An  unbiased  estimator  for  6^  is  the  crude  Monte  Carlo 
estimator  which  can  be  written  as 


(4.D 


Here  is  a  sample  from  the  cdf  Fx  which  is  generated  by  the  total  hazard 

construction  (3.3)-(3.5). 


Another  unbiased  estimator  for  0^  is 


(4.2)  [s(^W))+S(^"°)]  (k  is  even) 

where  X^m^  is  generated  by  (3.3)-(3.5)  using  the  vector  of  independent  uniform  [0,1] 
random  variables,  and  X^m^  is  generated  by  (3.3)-(3.5)  using  1  —  t/m\  We  call  5^  the 
antithetic  estimator  of  O^K 


4.  Some  applications  and  examples 


Since  g(X)  +  g(Xa)  is  a  special  case  of  g(X)  —  h(Y)  [with  h(Y)  replaced  by  —g(Xa)]  we  obtain 
from  Theorem  3.1  that  <  Far|0^].  That  is,  the  antithetic  estimator  is  more 

accurate  than  the  crude  Monte  Carlo  estimator 

Assume  now  that  we  seek  to  estimate  0®  =  J E[g(X)  -  h{Y)],  X,Y  e  R+,  i.e.,  the  expected  value 
of  response  difference  of  a  pair  of  functions  g(X)  and  h(Y).  Assume  that  both  X  and  Y  have 
CDTH  and  that  both  g  and  h  are  monotone  functions.  The  crude  Monte  Carlo  estimator  of 
0W  is 


(4.3)  0<2)  -  I J  [g(^m))  -  /KI(m))] 

m—  1 

where  X^m^  and  Y^m^  are  generated  independently  by  the  transformation  described  by  (3.3)-(3.5). 
The  antithetic  estimator  for  0(2>  is 

(4.4)  ef  - 1 2  -  Mlf0)] 

m— 1 

where  X  and  Y  are  generated  using  VCRN,  i.e.,  X^m^  is  generated  by  (3.3)-(3.5)  using 
and  Y^  is  generated  similarly  also  using  l/m\ 

It  follows  from  Theorem  3.1  that  <  Kar[0^]  i.e.,  the  antithetic  estimator  is  more 

accurate  than  the  crude  Monte  Carlo  estimator. 

The  monotonicity  assumption  in  Theorem  3.1  actually  holds  for  a  broad  class  of  complex 
stochastic  systems  (see,  e.g.,  Rubinstein,  Samorodnitsky  and  Shaked  (1985)).  Here  we  take 

g(x,,...,x5)  =  maxfminfxj,^),  min(xj,x3,^),min(x4,x5)), 
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h(ylt...,y5)  =  max(minO'1,^2),minCy1,^3,>'5),minO'4,>'5)). 


These  represent  the  lives  of  coherent  systems  described  by  a  Wheatstone  Bridge. 

Suppose  X=  has  hazard  rates  which  increase  with  the  number  of  failed  components. 

For  example,  for  a  >  0,  £  >  0,  let 

(4.5)  A,(x  |  &~  *,.•)  =  «(  Ul  +  fl./e/,/c  {1.....5}, 

where  \J\  denotes  the  number  of  elements  of  J.  Then  the  A/s  satisfy  (3.7)  and  thus  A- has 
CDTH. 

Let  y-(y,,...,ys)  have  hazard  rates  of  the  form  (4.5),  that  is,  for  f}>0,-q>0,  let  (using 
obvious  notation) 

(4-6)  Pj(y\  Yj  =  £j,‘)  =j3(|/|  +  T]),i  e.  J,J  c  £1, ...» 5}  . 

Then  also  Y  has  CDTH. 

For  the  choice  of  the  cdf  which  corresponds  to  (4.5)  it  is  easy  to  see  that 
A,  (x)  =  a£x,  x  >  0,  i  =  1 . 5, 

A|*  (x  |  Xj  =  Xj)  =  qj(  1 J  |  +  £)x,  x  >  0,/ e  J,/c  [1,...,5}  . 

Thus, 

Ar1(/)  =  K]"1r,  />  0,,=  1,..„ 5, 
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AnJ(/l^)  =  [a(UI  +  £)]  \  t>0,ieJ,Jc  ,5}  . 


These  were  used  in  (3.4)-(3.5)  to  generate  The  vectors  Y®  were  generated  similarly. 

The  parameters  0(1)  =  E[g{X )]  and  0(2)  =  E[g{X)  -  h{Y )]  were  estimated  for  various  choices 
of  a, /J,£  and  77  using,  respectively,  ^,6^,0^  and  8 £2).  The  estimators  8  and  8a  were 
calculated  according  to  (4.1)  and  (4.2)  while  estimating  8^  and  according  to  (4.3)  and  (4.4) 
while  estimating  0^2\ 

A  sensible  performance  measure  is  (/gVar(fl)  /  tg  War (8 a))  where  tg  and  tg  are  the  CPU  times 
for  the  crude  Monte  Carlo  and  the  antithetic  estimators  respectively.  Comparing  (4.1)  and 
(4.2)  [or  (4.3)  and  (4.4)]  we  see  that  tg  <  tg  because  the  antithetic  estimator  0a  needs  only 
half  as  many  random  numbers  as  it  is  needed  by  its  counterpart  6.  In  the  following  we  shall 
neglect  this  advantage  of  8a  and  take  the  performance  measure  to  be 

(4.7)  e  .  Var(0)  f  Var(0J. 
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Table  1 


The  efficiency  e  of  the  antithetic  estimator 
e?  as  a  function  of  a  and  £ 


1 

2 

3 

4 

5 

1 

2.05 

2.83 

2.63 

2.59 

2.27 

2 

1.90 

2.45 

1.88 

2.31 

2.40 

3 

2.21 

2.41 

2.16 

1.67 

1.58 

4 

2.10 

1.83 

1.93 

1.65 

2.05 

5 

1.81 

1.71 

1.79 

1.53 

2.18 

Table  1  gives  the  values  of  e  while  estimating  0(1*  for  various  choices  of  a  and  £.  Here,  in 
order  to  compute  0(1)  we  used  k  =  1000  replications  of  X.  For  the  computation  of  0^  we 
used  k/2  =  500  replications  of  the  pair  ( X ,  Xa ).  Since 


Kar(0(1))  =  Var(^g(X(Zm-1'>)  +  g(*(2/7,))). 


from  (4.7)  one  sees  that  e  =  Var(g{X +  g(X^2^))  /  Var(g(X('^)  +g(A^1)))-  Hence  e  is  esti¬ 
mated  by  e=*sz/s%,  where  s2  is  the  sample  variance  of  the  500  observations 

(g(^2m-1))+s(*<2m))), 

m  =  1,2, ...,  500  and  s%  is  the  sample  variance  of  the  500  observations 
(g(^m))  +  g(Xim])),  m  =  1 , 2, ....  500. 
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Table  2 


The  efficiency  of  e  of  the  antithetic  estimator 

em 

as  a  function  of  a,£,j 3  and  t] 


'Xtf.a) 

07.4)\ 

(2,1) 

(3,2) 

(4,3) 

(5,4) 

(6,5) 

(1,1) 

4.71 

14.73 

22.80 

38.63 

69.94 

(2,2) 

4.91 

15.08 

21.00 

37.28 

62.86 

(3,3) 

4.98 

14.88 

23.93 

39.53 

62.47 

(4,4) 

5.13 

11.51 

25.00 

43.01 

66.10 

(5,5) 

4.89 

11.00 

19.83 

36.20 

58.87 

Table  2  gives  the  values  of  e  while  estimating  6®  with  various  choices  of  a,£,j3  and  tj.  For 
-(2) 

the  computation  of  6  we  used  k  =  500  replications  of  (X,Y),  where  X  and  Y  are  independent. 

_(2) 

For  the  computation  of  0a  we  used  k  =  500  replications  of  ( X,  Ya).  Similar  to  the  earlier 
case  e  is  estimated  by  e  =  s2  /  52,  where  j2  is  the  sample  variance  of  the  500  observations 
g{X*m^)  —  h(Y^m^),m  =  1,2, ...,500  and  j2  is  the  sample  variance  of  the  500  observations 

g(X*m^)  —  m  =  1,2,..,  500. 

These  numerical  results  clearly  indicate  the  benefits  of  using  antithetic  and  common  random 
numbers  for  variance  reduction. 
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